It has been nearly 20 years since the first images of surface nanobubbles were reported by two independent groups working in China^[1] and Japan.^[2] Actually, the first clue to the existence of such surface nanobubbles can be traced from the earlier surface force measurements between hydrophobic surfaces immersed in aqueous solutions.^[3–5] As early as 1994, Parker, Claesson, and Attard proposed that the stepwise features in force curves from their experiments were due to the bridging of nanoscale bubbles adhered on solid surfaces.^[3] Following from these initial reports of surface nanobubbles, many experimental studies employing different techniques such as atomic force microscopy (AFM),^[1,2,6–15] rapid cryofixation,^[16] neutron reflectometry,^[17] and direct optical visualization^[18,19] have shown that surface nanobubbles can exist for a substantially long period of time.

Surface nanobubbles are of great interest as they have potential applications in many fields, such as boundary slip in fluid,^[20,21] froth-flotation,^[22] and protein adsorption.^[23,24] However, the stability of surface nanobubbles posed a theoretical conundrum.^[25–27] Theoretically, bubbles at the nanoscale are unstable owing to the sharp increase in the Laplace pressure. The Laplace pressure, Δp, for a spherical bubble of radius R is given by Δp = 2γ/R, where γ is the surface tension of the bubble interface. Therefore, smaller bubbles have higher internal pressure. For example, the extra internal pressure is ~ 14.3 atmospheres for a bubble of radius 100 nm in water. The increased internal pressure leads to an increase in the solubility of the gas in the solution surrounding the bubble. Consequently, if the solution around the bubble is not sufficiently oversaturated, gas will leave the bubble by diffusion in order to establish equilibrium. The loss of gas from the bubble leads to a decrease in the bubble size, which further raises the Laplace pressure and thus results in more gas leaving the bubble. As a result, a positive feedback cycle is established and rapidly leads to the dissolution and disappearance of small bubbles. The Epstein and Plesset theory,^[28] which was developed to describe the gas diffusion process around a bubble, predicts that in saturated solutions, small bubbles should very rapidly shrink and disappear. For example, the lifetime of bubbles smaller than 1000 nm is less than 0.02 s, making it is difficult to detect and measure such bubbles in experiments.^[29] Note that if the solution is oversaturated with gas, the positive feedback cycle will operate in reverse and will result in the rapid growth of bubbles from small bubbles, which also has been discussed by Epstein and Plesset.^[28] In short, the earlier studies based on classical theories predicted that there is no thermodynamic stability for surface nanobubbles, and such bubbles should shrink and disappear or grow to a macroscopic size in a very short time.

Different mechanisms, including the dynamic equilibrium theory^[30,31] and the contamination theory,^[32] have been proposed and provided important insights, but both of them are insufficient to account for stable surface nanobubbles.^[33] In 2012, we proposed that contact line pinning induced by surface heterogeneities (e.g., physical roughness or chemical heterogeneities) leads to the appearance of thermodynamically metastable surface nanobubbles, which could explain the existence of long-lived surface nanobubbles in experiments.^[34] The contribution of contact line pinning to the nanobubble stability has also been proposed in the theoretical work by Weijs and Lohse,^[35] and in the experimental work by Zhang et al.^[36] Several recent experimental works on surface nanobubbles also supported that contact line pinning is a necessary condition for achieving the stability of surface nanobubbles.^[37–41]

An impressive amount of experimental, theoretical, and computational work on surface nanobubble has been accumulated over the last 17 years. Here, we do not attempt to discuss the entire corpus of this field, for which the readers are encouraged to consult the more comprehensive reviews.^{[25–27,29,33,42–44]} Rather, we focus on recent theoretical and computational studies, especially those from our group, of pinned surface nanobubbles. From Section 2 to 5, we review related research on the pinned surface nanobubbles with different theories and methods, including the classical nucleation theory (CNT), Epstein and Plesset theory, lattice density functional theory (LDFT), and molecular dynamics (MD) simulation method. Finally, in Section 6, we draw conclusions and give a perspective on the present and future challenges.

On the basis of the pioneering work of Gibbs, Becker, and Doring,^[45] Volmer and Weber^[46] and others developed the CNT. According to the CNT, the free energy cost for the formation of a spherical bubble in bulk liquid consists of a volume term and a surface term: ΔG = −(4π/3)R³Δp + 4πR²γ, where Δp = p_v − p_l is the pressure difference between the center of bubble and the bulk liquid, γ is the surface tension, and R is the radius of the bubble. By setting ∂ΔG/∂R = 0, the free energy reaches its maximum at a critical radius R* = 2γ/Δp, which has the same formula as the Laplace equation for the mechanical stability of a bubble. Thus, bubbles will grow infinitely if their size exceeds the critical size, or they will disappear if they have a size less than the critical size.

When there exist surface heterogeneities that can pin the contact line of surface nanobubbles, the situation changes. In our previous work,^[34] we found that within the framework of CNT, a stable surface nanobubble corresponds to a free energy minimum when the contact line is pinned. The same conclusion was obtained by the subsequent theoretical work by Attard,^[47] in which the author used the modified CNT to include the oversaturation dependence of the surface tension.

CNT suggests that for stable surface nanobubbles, a relationship of sinθ = r/R, where θ is the contact angle of the nanobubble, r is the footprint radius of the bubble, and R is the surface radius that equals to the critical radius (i.e., R = R*), holds. The critical radius R* can be further written as R* = 2γ/Δp = 2γ/(ρ_vΔμ), where ρ_v is the density of vapor and Δμ is the difference in chemical potential between liquid and vapor, which determines the level of oversaturation. Therefore, this relationship predicts how the contact angle depends on the oversaturation and the footprint radius, as well as the upper threshold of the footprint radius (i.e., r < R*).

Within the framework of CNT, a negative feedback mechanism was proposed to explain the stability of pinned surface nanobubbles.^[34,48] The Laplace equation for a bubble in equilibrium with its surrounding, Δp = 2γ/R, describes the mechanical balance on the interface between the contribution from Δp, which makes the bubble grow, and that from γ, which makes the bubble shrink. Consider a pinned surface nanobubble as shown in Fig. 1(b), and a small perturbation of the bubble’s volume. If the volume of the nanobubble increases, the curvature radius of the bubble would decrease because the contact line cannot move, and therefore the contribution of 2γ/R increases. In this case, Δp < 2γ/R, such that the bubble shrinks accordingly and returns to the equilibrium state. If, conversely, the nanobubble shrinks in its volume under a small perturbation, its curvature radius would increase, and thus Δp > 2γ/R, leading to the growth of the bubble. Therefore, a negative feedback cycle is established to prevent the nanobubble from shrinking and growing, and stabilizes the nanobubble. This negative feedback mechanism only occurs when the contact angle of nanobubbles (measured in the liquid phase) is larger than 90°, which can explain why surface nanobubbles always have an increased contact angle compared with the corresponding angle of a macroscopic bubble or droplet on the same surface.^{[2,8,26,33,49]} However, if there is no contact line pinning, the equilibrium bubble is not stable and will rapidly disappear or grow under the positive feedback mechanism (see Fig. 1(a)).

Fig. 1.

	Figure Option ViewDownloadNew Window
	Fig. 1. (color online) Schematic illustration of the positive feedback mechanism for unstable nanobubbles without contact line pinning (a) and that of the negative feedback mechanism for stable nanobubbles with contact line pinning (b). The figure is reproduced from Ref. [48].

The kinetics of the gas diffusion process around a bubble were dealt in 1950 by Epstein and Plesset.^[28] A gas bubble in a solution will grow or shrink by diffusion, depending on the solution being oversaturated or undersaturated. To solve the diffusion equation, the local concentration of gas dissolved in the solution adjacent to bubble surface is determined by applying Henry’s law, by which the total pressure is the sum of the external pressure and the Laplace pressure. Therefore, the evolution of the bubble can be determined by the diffusion equation, Laplace equation, and Henry’s law.

For oversaturation (the gas oversaturation ξ > 0), Lohse and Zhang^[50] provided an exact calculation for a single pinned surface nanobubble and confirmed that the bubble can reach a stable equilibrium state, for which the bubble satisfies the relationship sinθ = ξL/L_c, where L is the bubble’s footprint diameter and L_c is the critical footprint diameter. This relationship can transform into that obtained from the CNT (see Section 2) by using L = 2r, L_c = 4γ/p₀(p₀ is the ambient pressure), Δp = ξp₀ (Henry’s law), and R* = 2γ/Δp (the critical radius in CNT). Therefore, above relationship also predicts the contact angle according to the oversaturation and the footprint diameter, as well as the upper threshold of the footprint diameter (i.e., L < L_c/ξ).

For undersaturation (ξ < 0), the theoretical work of Weijs and Lohse^[35] confirmed that pinned surface nanobubbles cannot reach a stable equilibrium state, but dissolve on a much longer timescale than free bubbles, which can explain why under normal experimental conditions (e.g., the liquid is exposed to atmospheric conditions), surface nanobubbles can live for many hours or even up to days.

The lattice density functional theory (LDFT),^[51–57] based on a simple cubic lattice gas model, provides a very simple but efficient method to investigate the behavior of simple fluids at a molecular level. LDFT has been widely used to study capillary condensation and evaporation,^[51–55] vapor–liquid nucleation,^[58–62] water bridges,^[63] and the wetting of solid surfaces.^[64,65] In recent years, LDFT has also been also used to study pinned surface nanobubbles.^[34,66–70]

4.1. Thermodynamically metastable pinned surface nanobubbles

LDFT was first used to explain the stability of surface nanobubbles in 2012.^[34] In that work, the free energy changes as a function of a nanobubble’s volume on various substrates were calculated by using the constrained LDFT.^[56] It was found that for oversaturated liquid on substrates with sufficient surface heterogeneities, surface nanobubbles can be in a thermodynamically metastable state: a state in which a nanobubble corresponds to a local minimum of free energy. Kinetic LDFT^[52] calculations showed that if the pinning effect were absent, the surface nanobubble would dissolve with roughly constant contact angle. However, if there existed the pinning effect, both the volume and the contact angle of the surface nanobubble would remain constant. Examples of stable pinned surface nanobubbles in oversaturated liquids, obtained by LDFT, are shown in Fig. 2.

Fig. 2.

Figure Option ViewDownloadNew Window

Fig. 2. (color online) Stable surface nanobubble induced by the contact line pinning effect. (a) From left to right, substrate with a ring pattern and stable surface nanobubbles at two different degrees of supersaturation. (b) Schematic drawing of the system. (c) Substrate with chemical heterogeneity in a ring-like shape and a stable surface nanobubble. (d) Solid surface with randomly distributed physical heterogeneity and stable surface nanobubble. The figure is reproduced from Ref. [34].

LDFT calculations also revealed that stable pinned surface nanobubbles always have a contact angle (measured from the liquid side) larger than 90°. Importantly, the contact angle was found to depend on both the oversaturation and the footprint radius of nanobubbles, but is independent of the substrate wettability.^[34] These results agree well with the theoretical results (see Sections 2 and 3) and experimental observations that the contact angle of nanobubbles is size-dependent and substrate-wettability-independent.^[2,8,26,33]

4.2. Relationship between substrates and pinned surface nanobubbles

LDFT was then used to investigate how the nanobubble stability depends on the substrate structure (roughness) and wettability.^[66] In the calculations, four types of substrates with different local structures (see Fig. 3) were considered: Type I (circular terrace), Type II (circular concave), Type III (ring), and Type IV (perforated ring) substrates. The substrate wettability was controlled by adjusting the strength of the fluid–solid interaction, ε_sf. Figure 3 also shows typical surface nanobubbles stabilized on those substrates at given oversaturation and pinning radius (i.e., the bubble’s footprint radius). Therefore, both the substrate roughness and wettability affect the nanobubble stability. However, these stable surface nanobubbles show the same morphology, indicating that the contact angle of a pinned surface nanobubble is independent of the substrate wettability, which again agrees with the theoretical predictions (see Sections 2 and 3).

Fig. 3.

	Figure Option ViewDownloadNew Window
	Fig. 3. (color online) Snapshots of four types of substrates (the first column) and stable surface nanobubbles (the second to fourth columns) on the substrates. From left to right, the substrate’s wettability is changed by adjusting the fluid–solid interaction strength, εsf. The figure is reproduced from Ref. [66].

To explain these results from LDFT calculations, the pinning force, f_p, which acts on the contact line of a nanobubble and prevents its lateral motion, was quantitatively evaluated according to the corresponding definition: f_p = γ(cosθ_bubble − cosθ₀), where γ is the surface tension, θ_bubble is the contact angle of the nanobubble, and θ₀ is the contact angle corresponding to a bubble on smooth a substrate.^[66] Both the pinning force required to stabilize nanobubbles and that which the substrate can provide were calculated. It was found that the substrate structure and wettability together determine the sign and threshold of the pinning force provided by the substrate, whereas the pinning force required to stabilize the nanobubble is related to the substrate wettability and the pinning radius. If and only if the required pinning force is within the range of the provided pinning forces, contact line pinning occurs and therefore the nanobubble is stable. Otherwise, the motion of the contact line occurs and the nanobubble becomes unstable.

Those calculations revealed the relationship between nanobubbles and substrates: “The substrate structure and wettability determine the possible range of pinning force and thus affect the nanobubble stability. However, for a stable nanobubble, the contact angle is independent of the substrate wettability.”

4.3. Interaction between AFM tips and pinned surface nanobubbles

In most experiments, the morphology of surface nanobubbles was imaged with AFM.^[1,2,6–15] The interaction between AFM tips and the surface nanobubbles certainly influences the final images of such bubbles. LDFT was also adopted to explore the process of the AFM tip approaching, contacting, and finally retracting from the pinned surface nanobubbles, and to unravel the interaction between them, as well as how the wettability and shape of the AFM tip affect the morphology of surface nanobubbles.^[69]

In the calculations, the system contained an AFM tip and a stable pinned surface nanobubble on a rough substrate. The rough substrate was modeled as a solid surface decorated by a ring-shaped pattern with a radius of 20σ (where σ is the lattice spacing), and the AFM tip was modeled as a hemispherical end on top of a solid cylinder with a radius of 10σ. A cone-shaped tip was also built to investigate the effect of tip shape.

LDFT calculations^[69] indicated that the nanobubble showed an elastic deformation in response to the approach of a hydrophilic tip (Fig. 4(a)). The hydrophilic nature of the tip and a thin wetting film covering it prevented the tip from penetrating the bubble during the approach process, showing an elastic effect. The elastic effect induces a strong nanobubble deformation and repulsive interaction. In the retraction process, however, the hydrophilic tip could be immersed in the bulk solution easily through departure from the nanobubble surface, leading to a weak and intermediate-range attraction between the tip and the nanobubble. For a hydrophobic tip, however, a different situation took place (see Fig. 4(b)). In this case, the nanobubble showed an adhering effect, which dominates the tip–nanobubble interaction. In the approach process, the hydrophobic nature of the tip facilitated its penetration through the vapor–liquid interface without causing significant nanobubble deformation, which results in the disappearance of the repulsive force. In the retraction process, the hydrophobic adhesion between the tip and the nanobubble induced a strong lengthening effect on the nanobubble deformation. As a result, a strong and long-ranged attractive force was observed.

Fig. 4.

	Figure Option ViewDownloadNew Window
	Fig. 4. (color online) Capillary force between AFM tip and nanobubble and the morphology of the nanobubbles for different tip–substrate distances in the approach (red arrow) and retraction process (black arrow): (a) for a hydrophilic tip and (b) for a hydrophobic tip. The figure is reproduced from Ref. [69].

However, LDFT calculations showed that a cone-shaped tip also induced the deformation of the nanobubble during the approaching and retracting process, though in a much weaker manner. Therefore, a sharp-ended hydrophilic tip maybe a good choice to design minimally invasive experiments, which was in good agreement with the conclusion of Walczyk and Schoönherr’s experimental work.^[71]

Compared to LDFT calculations, MD simulations are increasingly used in investigating surface nanobubbles.^[72–76]

5.1. Pinning and oversaturation induce stable surface nanobubbles

Stable pinned surface nanobubbles were first achieved by MD simulation in 2014.^[72] The stability of surface nanobubbles in both pure fluids and gas–liquid mixtures was studied. In the computational work, a fluid system with Lennard–Jones (LJ) potentials was employed to describe the interaction between the different species. The LAMMPS package^[77] was used, with constant temperature, pressure, and particle number. Note that in the MD simulations, the pressure was controlled by applying an external force on the wall at the top of the simulation box. A nanopore in the bottom substrate was explicitly included in the substrate to induce the contact line pinning. In the case with gas–liquid mixtures, a reservoir of gas molecules, in which the identity exchange of liquid and gas molecules in the reservoir was performed every 0.1 ns, was introduced to maintain a target gas concentration far from the surface nanobubble.

In the MD simulations,^[72] stable pinned surface nanobubbles were found in both pure liquids and gas–liquid mixtures, provided that there was suitable overheating or gas oversaturation (Figs. 5(a) and 5(c)). A Wenzel or a Cassie wetting state was found in undersaturated and saturated fluids, and the liquid-to-vapor phase transition occurred at high oversaturation (Figs. 5(b) and 5(d)). Moreover, the MD simulations proved that an excess of dissolved gas molecules is required to produce stable gas nanobubbles, which can explain why no nanobubbles are experimentally observed in degassed water at room temperature. However, for surface nanobubbles either in pure liquids or in gas–liquid mixtures, both the curvature radius and the contact angle of nanobubbles increase as the level of overheating/oversaturation decreases, which again agrees with the theoretical predictions (see Sections 2 and 3). These MD simulations of pinned surface nanobubbles^[72] unambiguously showed that the contact line pinning effect, together with the oversaturation turned out to be crucial for surface nanobubble stabilization.

Fig. 5.

Figure Option ViewDownloadNew Window

Fig. 5. (color online) (a) Stable pinned surface nanobubble in pure liquids at different degrees of overheating. From left to right, the degree of overheating is increasing. (b) Several typical snapshots during the liquid-to-vapor phase transition in pure liquids process at a high degree of overheating. (c) Stable pinned surface nanobubble in gas–liquid mixtures at different degrees of oversaturation. From left to right, the degree of oversaturation is increasing. (d) Several typical snapshots during the liquid-to-vapor phase transition process in gas–liquid mixtures at a high degree of oversaturation. The figure is reproduced from Ref. [72].

5.2. Nucleation of pinned surface nanobubble on rough substrates

MD simulations were then used to study how nanobubbles form on rough hydrophobic substrates.^[78] The simulation system employed in that work was similar to that in Ref. [72], except that the rough substrates were modeled as solid surfaces decorated with several identical nanopillars. Three different numbers of nanopillars were used to study the nanobubble nucleation at substrates with high, moderate, and low degrees of surface roughness. The substrates were hydrophobic, and the gas–liquid mixture was set to be oversaturated by controlling the gas concentration in a gas reservoir.

First, long-time standard simulations were carried out to directly observe the kinetic pathways.^[78] A two-step nucleation route involving the formation of an intermediate state was found for the nanobubble formation: first, several gas cavities occur in the grooves, leading to a Wenzel-to-Cassie transition (e.g., (i)–(iv) in Fig. 6(b)); then, those gas cavities coalesce and form a stable pinned surface nanobubble, inducing a Cassie-to-nanobubble transition (e.g., (iv)–(v) in Fig. 6(b)). Additionally, the corresponding free energy changes were quantitatively evaluated^[78] by using constrained simulations combined with the thermodynamic integration approach. The free energy barriers for the two sequential transitions have opposing dependencies on the degree of surface roughness, such that surface nanobubbles are more likely to form on the surfaces with moderate roughness as each state transition only needs to overcome a moderate free energy barrier.

Fig. 6.

Figure Option ViewDownloadNew Window

Fig. 6. (color online) Several typical snapshots during the nanobubble formations at surfaces with high roughness (a), moderate roughness (b), and low roughness (c). Snapshots (i)–(v) in panel (a), (i)–(v) in panel (b), and (i) in panel (c) are obtained from the standard long-time MD simulations, and the others are obtained from the constrained MD simulations using the thermodynamic integration approach. Reprinted with permission from Ref. [72].

5.3. Solvent exchange producing pinned surface nanobubbles

In experiments, the oversaturation environment that is required for stable pinned surface nanobubbles is often produced through the solvent exchange procedure. For example, Lou et al.^[1] used water (poor solvent for air) to exchange the ethanol (good solvent for air) to produce surface nanobubbles. MD simulations were recently used to understand the molecular mechanisms behind the solvent exchange.^[76] The simulation system was similar to that in Ref. [72]. The solvent exchange was achieved as follows: initially, the box was full of a good solvent; then, during the simulations, the identity exchange between liquid and gas molecules as well as the solvent exchange was performed at the same time interval, to maintain a stable gas concentration and solvent composition in a source region that was far from the substrate.

MD simulations demonstrated a two-stage mechanism for forming nanobubbles through a process of solvent exchange.^[76] During the first stage of the process, an interface between two interchanging solvents was found, which moves toward the substrate gradually as the exchange process proceeds. Unexpectedly, driven by the solubility gradient of liquid composition across the moving solvent–solvent interface, a directed diffusion of gas molecules against gas concentration gradient was found. The forced diffusion against the gas density gradient prevents the gas molecules from washing away and more importantly, it produces an increasingly high local gas oversaturation as the interface approaches the substrate. At the second stage, the local high gas oversaturation nucleates nanobubbles either on the solid surfaces or in the bulk solution, depending on the hydrophobicity of the substrate (see Fig. 7).

Fig. 7.

	Figure Option ViewDownloadNew Window
	Fig. 7. (color online) (a)–(c) Typical snapshots at different simulation time and (d)–(f) time evolution of the local gas density in different regions during the solvent exchange processes. The Young’s contact angle for the substrates are 130° ((a), (d)), 91° ((b), (e)), and 31° ((c), (f)), respectively.

The origin of the stability of surface nanobubbles has been a controversial topic for a long time. Now, an increasing body of experimental, theoretical, and computational work has confirmed that contact line pinning, which results from the intrinsic nanoscale physical roughness or chemical heterogeneities of substrates, is crucial for the nanobubble stability.

Contact line pinning implies there exists the negative feedback mechanism to stabilize surface nanobubbles: for increasing bubble volume, the Laplace pressure that makes the bubble shrink increases accordingly, very different from free nanobubbles, for which the Laplace pressure decreases once the bubble volume increases. Under contact line pinning, both the classical nucleation theory and the classical Epstein and Plesset theory^[28] predicted that surface nanobubbles could reach a thermodynamically stable state in oversaturated liquids/solutions,^[34,50] following a simple size constraint, i.e., sinθ = r/R and R = R* = 2γ/Δp = 2γ/(ρ_vΔμ). This size constraint shows that the contact angle only depends on the oversaturation and the footprint radius of the nanobubble, and can explain the anomalous contact angle for surface nanobubbles observed in experiments. These conclusions were further confirmed by using the LDFT^[34] and MD simulations.^[72]

Following from these discussions for the stability of pinned surface nanobubbles, both LDFT and MD simulations were widely utilized to reveal their various behaviors, such as the pinning force,^[66] interaction between AFM tips and nanobubbles,^[69] and nanobubble formation.^[76,78]

Although those studies give us useful insights into pinned surface nanobubbles, several challenges and limitations in these calculations and simulations still persist. The lattice nature embedded in LDFT calculations limits its ability to provide finer molecular details. The single-component model employed in present LDFT calculations fails to account for the role of dissolved gas in surface nanobubbles. Although stable surface nanobubbles exist in both pure liquid and gas–liquid mixtures,^[72] dissolved gas was found to have a significant influence on the bubble growth dynamics.^[79] Thus, it would be of interest to use the LDFT extended to a binary mixture model^[54] in future research. The effect of dissolved gas can be easily included in MD simulations. However, owing to the limit of computational resources, in most present MD simulations of pinned surface nanobubbles, only LJ potentials have been used, with length scales of at most tens of nanometers, and time scales of at most hundreds of nanoseconds. It would also be of interest to use more real models and perform longer simulations to gain simulation data that could be directly compared with experimental observations. Of course, other computational methods should be considered for application in this field. For instance, the lattice Boltzmann method,^[80] which is an approach in between continuum dynamics simulations and MD simulations, may have great potential to study the collective effect^[81] in the dissolving and formation process of surface nanobubbles.^[33] This type of study would certainly deepen our understanding of cavitation and bubble nucleation and collapse.^[82–87]